Fast crack propagation by surface diffusion 
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We present a continuum theory which describes the fast growth of a crack by surface diffusion. 
This mechanism overcomes the usual cusp singularity by a self-consistent selection of the crack tip 
radius. It predicts the saturation of the steady state crack velocity appreciably below the Rayleigh 
speed and tip blunting. Furthermore, it includes the possibility of a tip splitting instability for high 
applied tensions. 
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One of the most challenging puzzles in nonequilibrium 
physics and materials science is the phenomenon of frac- 
ture. It is important for the vast field of material fail- 
ure and probably also for friction processes [^j. Despite 
its relevance even the motion of a single crack is poorly 
understood Q|. Experimentally, the maximum attained 
crack velocities are far lower than the theoretically ex- 
pected Rayleigh speed Beyond a critical velocity, a 
so far unpredictable tip splitting of the crack can happen 
and produce strange oscillations of the crack speed Q . 

The classical theories are based on linear theory of 
elasticity and on an integral energy balance in the vicin- 
ity of an infinitely sharp crack tip. However, a more 
detailed approach based on equations of motion for the 
crack shape is needed to describe the intriguing spectrum 
of phenomena near the crack tip. Hence, the curvature 
of the crack tip is required as a new relevant dynami- 
cal variable which also allows to avoid stress singulari- 
ties. We emphasize that, in contrast to models which 
describe crack propagation by bond breaking at the in- 
finitely sharp tip, growth with a finite tip radius always 
requires a transport mechanism in order to preserve the 
shape (see Fig. |l|). Recently it was proposed that this 
lengthscale is dynamically selected by the threshold of 
plastic deformations in the tip region Unfortunately, 
approaches of this type (see also |Q) require the intro- 
duction of theories of plasticity which are usually much 
more speculative and less verified than the ordinary lin- 
ear theory of elasticity. 

Here we demonstrate that linear theory of elasticity 
is sufficient to describe consistently crack propagation, 
driven by surface diffusion along the crack surfaces. Of 
course, in many situations plasticity is very important, 
but the beauty of our approach is that it predicts, in a 
simple and well controlled continuum theory, steady state 
crack growth, the tip splitting instability and also slow 
deformations of already existing cracks. 

The idea of crack propagation by surface diffusion has 
previously been studied by Stevens and Button j?), who 
assumed an ad hoc crack shape which is not found by 
solving the full free boundary problem and requires mass 
transport over large scales. Therefore, their model can- 
not describe the usual fast crack growth. 



Our basic mechanism is related to the Asaro-Tiller- 
Grinfeld (ATG) instability which predicts a morpho- 
logical instability of a uniaxially stressed solid interface 
due to surface diffusion. Relatively long-wave perturba- 
tions of the interface lead to a reduction of the elastic 
energy of the system, whereas short-wave corrugations 
are hampered by surface energy. In the long time be- 
havior, deep grooves can form, producing shapes similar 
to cracks According to previous theories, which 

used only the static theory of elasticity, the notches prop- 
agate with increasing velocity and decreasing tip radius 
and collapse to a finite-time cusp singularity. Similar to 
crack dynamics the lack of tip radius selection becomes 
obvious, and already shows the close relationship between 
the ATG instability and crack propagation. 

Usually, it is believed that surface diffusion is slow, 
but, surprisingly enough, it should not be ignored even 
in fast fracture processes. Our main idea is that surface 
diffusion is driven by the strong gradient of the chemi- 
cal potential in the tip region. This can be a very effi- 
cient mechanism for crack propagation if the transport 
length is sufficiently small. Additionally, energy release 
and strong dissipation bring the local temperature close 
to the melting temperature |11, |2|. This drastically in- 



creases the surface diffusion coefficient and makes fast 
crack propagation essentially independent of the outside 
temperature. 

Steady state crack growth — We use a two-dimensional 
plane-strain situation with mode I loading to describe 
crack propagation |^ . 

On the surfaces of the crack the normal stress cr„„ and 
the shear stress (t„t- vanish, whereas the tangential stress 
arr usually does not. The chemical potential at the in- 
terface is given by 



an 



(1) 



Here a is the surface energy, n the curvature of the in- 
terface and SI the atomic volume. E and v are Young's 
modulus and Poisson ratio, respectively. Non-hydrostatic 
stresses drive a surface flux proportional to the gradients 
of the chemical potential along the surface; in turn the 
normal velocity equals the divergence of this flux due to 
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conservation of material, 



D 



(2) 



where d/ds denotes the tangential derivative and D (di- 
mension m'*s^^) is proportional to the surface diffusion 
coefficient. (It is related to the usual surface diffusion 
coefficient Ds hy D = Dg^^Sa/kT. Here 5 is the num- 
ber of atoms per unit area of surface, k the Boltzmann 
constant and T the temperature.) 

First, we are interested in steady state solutions of the 
equation of motion, with a crack moving in positive x- 
direction with velocity v (see Fig. 1); In co-moving po- 
lar coordinates, x = r{9) cos 9, y — r{9) sinO, the steady 
state equation for the shape r{9) reads after one integra- 
tion of Eq. (|) 



vr sin 9 



D 



an Vr^ + r'2 d9 ' 



(3) 



Generally speaking, this, together with Eq. (|I|), is a com- 
plicated, non-linear third order equation with non-local 
contributions arising from the elastic fields, since ctxt de- 
pends on the entire shape. 

In the tail region stresses decay and the shape equa- 
tion is a third order linear differential equation Dy'" = vy 
with two growing (and oscillating) and one decaying so- 
lution. Only the latter, y{x —oo) — Aexp[{v/Dy^^x], 
asymptotically describes physical shapes and is allowed. 
Let us focus on symmetrical solutions, r{d) = r{—9), 
and start integration at the crack tip ^ = 0. Since the 
physical properties, curvature and stresses, do not de- 
pend on the choice of coordinate system but only on the 
crack shape, we can arbitrarily chose r{0 = 0) = rg, 
with the a priori unknown tip radius — 1/k(0). Then 
from symmetry and the definition of the tip curvature, 
K= {r'^ -\-2r''^ — rr")/{r'^ +r''^Y/'^, the natural conditions 
r'(0) = r"(0) = arise. Integration over the upper inter- 
face 9 > Q requires the suppression of two growing expo- 
nentials at the tail, which imposes two boundary condi- 
tions. For a given external loading, these conditions can 
be fulfilled by a proper selection of the tip radius and 
growth velocity v. By this argument the situation seems 
to be fully described. However, as we mentioned already 
earlier, the use of a static theory of elasticity does not 
allow a selection of the tip radius. The reason is, that 
both contributions to the chemical potential, surface en- 
ergy ~ K and elastic energy ~ cr^, behave as r^^ 
close to the tip; In the tip approximation, stresses behave 
as i 



(4) 



with the static stress intensity factor K ~ aocL^^"^, where 
tJoo is the applied remote stress and L is the macroscopic 



length of the crack, L ^ ro, which is not considered here. 
Instead we assume K to be kept fixed. The universal 
stress distribution depends only on the orientation 
relative to the crack The asymptotic distribution 
(^) is valid far away from the tip, tq <C r ^ L, but 
it gives the correct scaling of stresses also on the crack 
surface r « tq. Therefore, a dimensionless rescaling of 
all lengthscales, e.g. r = r/ro, and of the growth velocity 
V — vvq/D leaves the equation of motion invariant and, 
thus, cannot determine the lengthscale ro- Consequently, 
a steady state solution does not exist. This is the reason 
for the already mentioned cusp singularity of the ATG 
instability. 

The main idea of this paper is based on the fact that 
a full elastodynamic description restores the selection of 
this lengthscale. It is known that at least for higher crack 
speeds the angular distribution fij become strongly de- 
pendent on the ratio v/vr {vr is the Rayleigh speed [^). 
The dynamical stress intensity factor Kdyn is related to 
the static one used here by an extra velocity dependent 
function g{v/vB), Kdyn — Kg{v/vR). The crucial obser- 
vation is, that velocity appears now in two different com- 
binations in the equation of motion, vvq/D and v/vr. 
Thus, by introduction of the new parameter v/vr, a se- 
lection of both V and ro happens. 

From these general arguments we conclude that fast 
steady state crack propagation by surface diffusion is in- 
deed possible. However, the exact solution of the problem 
is technically very difficult, because it requires the solu- 
tion of an elastodynamic problem for an a priori unknown 
crack shape. The bulk equations of elasticity 



do 



dxj 



(5) 



are subject to the boundary conditions on the crack sur- 
face (surface of discontinuity) H 



-f piliVn = 



(6) 



with the mass density p and the elastic displacement field 
Ui. This is just the momentum balance equation on the 
free surface which moves with normal velocity u„. Fi- 
nally, the expression for the chemical potential on the 
crack surface should also be corrected compared to the 
static case: 



p = n[ ^(TikUik - ^piif - an 



(7) 



Eqs. (||)-(0) together with surface diffusion Eq. (||) and 
given loading configuration describe the crack propaga- 
tion in our model. We note that inertial effects appear 
not only in the bulk equations of elasticity but also in the 
boundary conditions and in the expression for the chem- 
ical potential. All these effects lead to the appearance of 
the parameter {v/vrY in the problem compared to the 
quasistatic description. 
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The preceding equations can be derived from 

the Lagrangian 



L = 



V(t) 



a. 

2 



dV - 



TtdS. 



S(t) 



(8) 



with V{t) being the time dependent volume of the sohd 
and S{t) its surface. The elastic equations and boundary 
conditions follow by the condition that S = J Cdt is sta- 
tionary with respect to variations of the displacement w^, 
and the chemical potential is related to the variation of 
C with respect to the interface position Evalu- 
ation of the exact equations of motion requires extended 
numerics, especially in the time-dependent case. 

The local crack tip model — We simplify the problem in 
order to make further analytical progress and to expose 
the general idea of our approach. It will turn out that 
one cannot describe all effects by this approximation and 
further refinement is necessary, but the main results are 
qualitatively very robust against changes of the model. 

We mimic the tangential stress by a local description in 
the spirit of Eq. (^) , as depending on the propagation ve- 
locity and only on the local properties of the interface. It 
takes both the velocity dependence of the angular distri- 
bution and the decrease of the dynamical stress intensity 
factor into account: 



K 



^l-{v/vRYcos{e/2) + {v/vRfsmH\ jv^l^. 

(9) 

This form reflects the first order transition of the princi- 
pal stress direction = for low velocities towards 7^ 
as function oivjvR |^. The use of more sophisticated ex- 
pressions (e.g. the singular dynamical field in full detail) 
would not provide a large gain, since, anyway local ap- 
proximations cannot lead to exact results. For the same 
reasons we also neglect inertial corrections to the bound- 
ary conditions in (H) and in the chemical potential (^. 

However, we have checked that a model with a contin- 
uous transition in azimuthal stress (replacement of sin* Q 
by sin^ Q in Eq. |^) gives qualitatively the same results. 

Now Eq. (|), together with ^ and (|), is a closed 
third order differential equation for the shape r(ff) which 
can be easily integrated numerically. It provides both 
the crack shape (Fig. ^ and a selection of vjvR and 
ro as functions of the dimensionless driving force A = 
K^{ \ — v'^)/2Ea. The results are given in Figures ^ and 

One of the main results is that the upper limit for 
the steady state crack velocity is appreciable below the 
Rayleigh speed, as known from experimental results. The 
instantaneous velocity in the non steady state regime can 
of course reach higher values For relatively low 

driving forces, the growth velocity increases with increas- 
ing A, but for higher values of A, it even decreases. Si- 
multaneously, the tip becomes sharper at first, but then 
blunts again. 
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FIG. 1: Calculated shape of the crack (without elastic dis- 
placements) driven by surface diffusion for A = 2. The ad- 
vance of the crack in positive x direction is indicated by the 
dashed curve. This requires the redistribution of matter along 
the crack by a transport mechanism. 




FIG. 2: Steady state velocity of the crack versus dimensionless 
driving force A. 
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FIG. 3: Dimensionless crack tip radius ro versus dimensionless 
driving force A. 
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Stability — Although the decrease of the velocity v as 
function of A might be naively understood as a sign of 
instability, the model itself is stable: We performed a 
straightforward but tedious numerical stability analysis 
and found no unstable modes. The point is, that only 
v/vr decreases, but vr^/D increases. Only a decrease of 
this parameter, which appears in combination with the 
dissipative coefficient D, would be a real sign of instabil- 
ity. 

Nevertheless, our solution is subject to the ATG insta- 
bility above a critical threshold of the driving force A. 
The de-stabilizing effects stem from the non-local elas- 
tic contributions which we neglected in our model, but 
which are obviously present in the real problem. This 
can be proved by simple dimensional analysis arguments: 
The characteristic wavelength of the ATG instability is 
A - Ea/al^ H; in the tip region it reads A^^p) ro/A. 
Thus, as soon as a certain critical driving force Ac is 
exceeded, the characteristic wavelength of instability fits 
into the tip region. The material independent number 
Ac is the threshold for the instability of the steady state 
solution, whereas A = 1 is the Grifhth point. Since, ac- 
cording to the steady state solution, v/vn is a universal 
function of the dimensionless parameter A, the threshold 
of instability in terms of v /vr is also essentially material 
independent. It is important to note that our steady 
state predictions are valid only below the threshold of 
instability. Thus, the main part or all of the decrease of 
velocity versus A is screened by the instability. Beyond 
the instability point, the behavior of the system is gov- 
erned by the full time dependent evolution. Hence we 
expect the ATG instability to be the relevant mechanism 
for the experimentally observed microbranching instabil- 
ity iQ, In contrast to the long wave instability [ p^ , 
this instability is localized in the tip region and cannot 
be suppressed by convective effects. 

Conclusion — We have developed a self-consistent con- 
tinuum model for crack propagation in homogeneous me- 
dia. Both ingredients of our theory, the linear elasticity 
which is valid everywhere in the bulk and surface dif- 
fusion which provides a mass transport and dissipative 
mechanism for crack propagation, are well established. 
The model is essentially parameter free, leading to the 
prediction that two dimensionless quantities, the crack 
velocity - v/vr and crack tip radius - ro(vii/Dy^^ are 
universal functions of the dimensionless driving force A. 
Strictly speaking, these functions still depend weakly on 
the Poisson ratio. We note that these statements, to- 
gether with the prediction of the tip instability above 
some critical velocity, are based on the general structure 
of our theory and do not involve the specific modeling 
of the surface stresses. The specific results given in the 
figures should differ from exact solutions only quantita- 
tively. 

It is important to realize that our model does not con- 



tradict classical theories g], but contains more informa- 
tion. This allows to calculate both the crack velocity and 
the fracture energy, while the classical theories predict 
only a relation between these quantities (the integral en- 
ergy balance). In our model, this energy balance is not 
fulfilled by approaching the Rayleigh speed limit with 
increasing driving force, as the classical theories predict, 
but by an increase of the fracture energy by tip blunting, 
which eventually leads to the tip instability. 

The scale of velocity is set by vr and is independent 
on the diffusion coefficient. However, the length scale of 
the tip is set by {D /vr)'^/^ . As we already noted, strong 
dissipation brings the local crack tip temperature close 
to the melting temperature. The diffusion coefficient Ds 
is then about 10~^m^s~^ and independent of outside 
temperature. Consequently, the lengthscale of the tip is 
of order of atomic units. We hope, that even for these 
small scales our general, qualitative predictions remain 
correct. 
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